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Abstract 

We evaluate the relic density of neutralinos in the minimal supergravity (mSUGRA) model. All 2^2 
neutralino annihilation diagrams, as well as all initial states involving sleptons, charginos, neutralinos and 
third generation squarks are included. Relativistic thermal averaging of the velocity times cross sections 
is performed. We find that co-annihilation effects are only important on the edges of the model parameter 
space, where some amount of fine-tuning is necessary to obtain a reasonable relic density. Alternatively, 
at high tan/3, annihilation through the broad Higgs resonances gives rise to an acceptable neutralino relic 
density over broad regions of parameter space where little or no fine-tuning is needed. 


1 Introduction 


A wide variety of astrophysical measurements are being used to pin down some of the basic cosmological 
parameters of the universe. High resolution maps of the cosmic microwave background (CMB) radiation|Q 
imply that the energy density of the universe Ul = pf Pc — 1, consistent with inflationary cosmology. 
Here, pc = 3F[^/8ttGn is the critical closure density of the universe, where Gn is Newton’s constant 
and H = lOO/i km/sec/Mpc is the scaled Hubble constant. The value of h itself is determined to be 
^ 0.7±0.1 by improved measurements of distant galaxies^. Meanwhile, data from distant supernovae|§ 
imply a nonzero dark energy content of the universe Ha ~ 0.7, a result which is confirmed by fits 
to the CMB power spectrum[Q. Analyses of Big Bang nucleosynthesis!^ imply the baryonic density 
Tlbh^ ~ 0.020 ± 0.002, although the CMB fits suggest a somewhat higher value of ~ 0.03. Hot dark 
matter, for instance from massive neutrinos, should give only a small contribution to the total matter 
density of the universe. In contrast, a variety of data ranging from galactic rotation curves to large scale 
structure and the CMB imply a significant density of cold dark matter (CDM)j|] ~ 0.2 ± 0.1. 

In many i?-parity conserving supersymmetric models of particle physics, the lightest neutralino (Zi) 
is also the lightest SUSY particle (ESP); as such, it is massive, neutral and stable. For this case, 
relic neutralinos left over from the Big Bang provide an excellent candidate for the CDM content of 
the universe [Q. In this work, we present results of calculations of the neutralino relic density within 
the context of the paradigm minimal supergravity model (mSUCRA, or CMSSM)||]. In mSUCRA, it 
is assumed that SUSY breaking occurs in a hidden sector of the model, with SUSY breaking effects 
communicated from hidden to observable sectors via gravitational interactions. The model parameter 
space is given by 

itt-Oj 'm-i/ 2 , Aq, tan/3 and sign{p). (1) 


Here, mg is the universal scalar mass, mi /2 is the universal gaugino mass and Aq is the universal trilinear 
mass all evaluated at Mgut, while tan /3 is the ratio of Higgs field vevs (vu/vd), and ^ is a supersymmetric 
Higgs mass term. The soft SUSY breaking parameters, along with gauge and Yukawa couplings, evolve 
from Mgut to Myjeak according to their renormalization group (RC) equations. At M^jeak, the RC 
improved I-loop effective potential is minimized, and electroweak gauge symmetry is broken radiatively. 
In this report, we implement the mSUCRA solution encoded in ISAJET v7.64[^. 

There is a long history of increasingly sophisticated solutions for the relic density of neutralinos in 
supersymmetric models [|l0|, 0 H, H, S 0 0,0,0 0 0,0,0 0 0, 


H 


35 . The key ingredient to solving the Boltzmann equation is to evaluate the thermally 
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averaged neutralino annihilation cross section times velocity factor. Traditionally, the solution is made 
by expanding the annihilation cross section as a power series in neutralino velocity, so that angular and 
energy integrals can be evaluated analytically. The remaining integral over temperature can then be 
performed numerically. The power series solution is valid in many regions of model parameter space 
because the relic neutralino velocity is expected to be highly non-relativistic. 

However, it was emphasized by Griest and Seckel that annihilations may occur through s-channel reso¬ 
nances at high enough energies iQ that a relativistic treatment of thermal averaging might be necessary. 
Drees and Nojiri found that at large values of the parameter tan /3, neutralino annihilation can be dom¬ 
inated by s-channel scattering through broad A and H Higgs resonances ||T^. The proper formalism for 
relativistic thermal averaging was developed by Gondolo and Gelmini (GG)[^, and was implemented 
in the code of Baer and Brhlik|^ Working within the framework of the mSUGRA model, it was 
found m, m, Q that at large tan/3, indeed large new regions of model parameter space gave 
rise to reasonable values for the CDM relic density. At large tan /3, the A and H resonances are broad 
enough (typically 10-50 GeV) that even if the quantity 2m~ is several partial widths away from exact 
resonance, there can still be a significant rate for neutralino annihilation. Thus, in the mSUGRA model 
at low mg and tan /3, neutralino annihilation is dominated by t-channel slepton exchange, and reasonable 
values of the relic density occur only for relatively low values of mg and mi/ 2 - At high tan/3, a much 
larger parameter space is allowed, owing to off-resonance neutralino annihilation through the broad Higgs 
resonances. 


In addition, there exist regions of mSUGRA model parameter space where co-annihilation processes are 
important, and even dominant. It was stressed by Griest and Seckel [pT| that in regions with a higgsino- 


like LSP, the Zi, Wi and Z 2 masses become nearly degenerate, so that all three species can exist in 
thermal equilibrium, and annihilate against one another. The relativistic thermal averaging formalism of 
GG was extended to include co-annihilation processes by Edsjo and Gondolo^^, and was implemented 
in the DarkSUSY code|^ for co-annihilation of charginos and heavier neutralinos. 


The importance of neutralino-slepton co-annihilation was stressed by Ellis et al. and others ||26|, P7| , 
p0| . In regions of mSUGRA parameter space where Zi and fi (or other sleptons) were nearly 
degenerate (at low mg), co-annihilations could give rise to reasonable values of the relic density even at 
very large values of mi/ 2 , at both low and high tan/3. In addition, for large values of the parameter Ag 
or for non-universal scalar masses, top or bottom squark masses could become nearly degenerate with 
the Zi, so that squark co-annihilation processes can become important as well[^ 


The goal of this study is to calculate the relic density of neutralinos in the mSUGRA model including 
co-annihilation processes in addition to relativistic thermal averaging of the annihilation cross section 
times velocity. Since there are very many Feynman diagrams to evaluate for neutralino annihilations and 
co-annihilations, we use CompHEP v.33.23|Q, which provides for fast and efficient automatic evaluation 
of tree level processes in the SM or MSSM. For initial states including Zi, Z 2 , Wi, ei, p,i, fi, ii and bi, 
we count 1722 subprocesses, including 7618 Feynman diagrams. For those processes we have calculated 
the squared matrix element and have written it down in the form of CompHEP FORTRAN output. 


The weak scale parameters from supersymmetric models are generated using ISAJET v7.64, and inter¬ 
faced with the squared matrix elements from CompHEP. Details of our computational algorithm are given 
in Sec. 2. In Sec. 3 we present a variety of results for the relic density in mSUGRA model parameter 
space. Much of parameter space is ruled out at low tan /3 since the relic density is too high, and would 
yield too small an age of the universe. At high tan /3, large regions of parameter space are available with 
a reasonable relic density in the range 0.1 < D~ < 0.3. In Sec. 4, we conclude. 


2 Calculational Details 

The evolution of the number density of supersymmetric relics in the universe is described by the Boltz¬ 
mann equation as formulated for a Friedmann-Robertson-Walker universe. For calculations including 
many particle species, such as the case where co-annihilations are important, there is a Boltzmann equa¬ 
tion for each particle species. Following Griest and SeckellQ, the equations can be combined to obtain 
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a single equation 


dn 

dt 


N 

—3Hn — (aeffv) (n^ — , where n = m 

i=l 


( 2 ) 


and the sum extends over the N particle species contributing to the relic density, with rn being the 
number density of the zth species. Furthermore, neq,i is the number density of the *th species in thermal 
equilibrium, given by 


^eq,i 


27r^ 




(3) 


where Kj is a modihed Bessel function of the second kind of order j. 

The quantity {ueffv) is the thermally averaged cross section times velocity. A succinct expression 
for this quantity using relativistic thermal averaging was computed by Gondolo and Gelmini for the 
case of a single particle species |j^, and was extended by Edsjo and Gondolo for the case including 
CO- annihilations!^. We adopt this latter form, given by 


(CTe//u)(x) 


ir -^1 (f) T,ij=i A(a^bf,6j)g^gJg,J-(a)da 


(4) 


where x = Tim— is the temperature in units of mass of the relic neutralino, ( 7 ^ is the cross section for 
the annihilation reaction ij ^ A (A is any allowed final state consisting of 2 SM and/or Higgs particles), 
X{a'^, bf, bj) = a'^ + bf + bj — 2(a^6f -b -b bfb'j), a = ^jm- and bi = mijm- . This expression is our 
master formula for the relativistically thermal averaged annihilation cross section times velocity. 

To solve the Boltzmann equation, we introduce a freeze-out temperature Tp, so that the relic density 
of neutralinos is given bjfl 

-4 (5) 

8.1 X 10-47 GeV"* 


where 


9(T.) 1.66 A 



rp3 




10“^ {^effv)dx 


( 6 ) 


The freeze-out temperature xp = Tp jm— is determined as usual by an iterative solution of the freeze-out 
relation 


X 


-1 

F 


log 


9eff 
27r3 2 


45 


2gtfGN 


{aeffv){xp) 



(7) 


Here, geff denotes the effective number of degrees of freedom of the co-annihilating particles, as defined 
by Griest and SeckelQ. The quantity 5 * is the SM effective degrees of freedom parameter with ~ 9 
over our region of interest. 

The challenge then is to evaluate all possible channels for neutralino annihilation to SM and/or Higgs 
particles, as well as all co-annihilation reactions. The 7618 Feynman diagrams are evaluated using 
GompHEP. To achieve our hnal result with relativistic thermal averaging, a three-dimensional integral 
must be performed over i.) the final state subprocess scattering angle 9, ii.) the subprocess energy 
parameter a = ^Jsjm— , and in.) the temperature T from freeze-out Tp to the present day temperature 
of the universe, which can effectively be taken to be 0. We perform the three-dimensional integral 
using the BASES algorithmic^, which implements sequentially improved sampling in multi-dimensional 
Monte Garlo integration, generally with good convergence properties. We note that the three-dimensional 
integration appearing in the case of our relativistic calculations involving several species in thermal 
equilibrium is about 2 orders of magnitude more GPU-time consuming than the series expansion approach, 
which requires just one numerical integration. 

®The procedure we follow gives numerical results valid to about 10% versus a direct numerical solution of the Boltzmann 
equation|M]. 
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3 Results 


Our first results in Fig. show regions of 0~ h?' in the toq vs. toi/2 plane in the minimal supergravity 
model for Aq = 0, tan/3 = 10 and for /r < O(left) and /i > O(right). 
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Figure 1: Regions of neutralino relic density in the mo vs. mi /2 plane for Ag = 0 and tan/3 = 10. 


The dark shaded (red) regions are excluded by 
theoretical constraints (lack of REWSB on the right, 
a charged LSP in the upper left). The unshaded ^ 
regions have > 1, and should be excluded, 

as they would lead to a universe of age less than -g 
10 billion years, in conflict with the oldest stars ^ 
found in globular clusters. The medium shaded ^ 
(green) region yields values of 0.1 < h? < 0.3, b 

. V 

i. e. in the most cosmologically favored region. The ^ ^ 
light shaded (yellow)< 0.1) and black(blue) ^ 
(0.3 < ^ correspond to regions with in¬ 

termediate values of low and high relic density, re¬ 
spectively. Points with mi /2 ~ 150 GeV give rise 
to chargino masses below bounds from LEP2; the 
LEP2 excluded regions due to chargino, slepton and 
Higgs searches are not shown on these plots. The 
structure of these plots can be understood by ex¬ 
amining the thermally averaged cross section times 
velocity, integrated from zero temperature to Tp. 

In Fig. 1^ we show this quantity for a variety of con¬ 
tributing subprocesses plotted versus mg for fixed 
Wi /2 = 300 GeV, A* > 0, and all other parameters as 
in Fig. At low values of mg, the neutralino annihi¬ 
lation cross section is dominated by t-channel scat¬ 
tering into leptons pairs, as shown by the black solid 
curve. However, at the very lowest values of mg, the 
annihilation rate is sharply increased by neutralino- 
stau and stau-stau co-annihilations, leading to very 



Figure 2: Thermally averaged cross section times 
velocity integrated from T = 0 to Tp, for var¬ 
ious subprocess. The thick light-grey(light-blue) 
curve denotes the total of all annihilation and co¬ 
annihilation reactions; mj ^/2 = 300 GeV, /i > 0, 
Ag = 0 and tan/3 = 10. 
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low relic densities where m~ ~ mf, As toq increases, the slepton masses also increase, which sup¬ 
presses the annihilation cross section, and the relic density rises to values 17— h? > 1. When mg increases 
further, to beyond the ~ 1 TeV level, and approaches the excluded region, the magnitude of the /r pa¬ 
rameter falls, and the higgsino component of Zi increases. This is the so called “focus point” region, 
explored in Ref. |^. In this region, the annihilation rate is dominated by scattering into WW, ZZ, 
hh and Zh channels. At even higher mg values, m~ ~ to— ~ to— , and these co-annihilation channels 
increase even more the annihilation rate. Finally, at the large TOg bound on parameter space, |/i| ^ 0, 
and appropriate REWSB no longer occurs. Most of the structure of Fig. || can be understood in these 
terms, with the exception being the horizontal band of very low relic density at TOi /2 — 125 GeV. In this 
region, which is nearly excluded by LEP2 bounds on the chargino mass, there is enhanced neutralino 
annihilation through the Z and h resonances. In fact, a higher degree of resolution on our plots would 
resolve these horizontal bands into two bands, corresponding to each of the separate resonances, as shown 
in Ref. |2^. The TOg vs. mi /2 planes for tan/3 = 30 are shown in Fig.|[ The structure of these plots are 
qualitatively the same as in Fig. |^. Quantitatively, they differ mainly in that the cosmologically favored 
regions are expanding as tan/3 grows. One reason is that t^e^light stau becomes even lighter as tan/3 
increases, and this increases the neutralino annihilation rate .^1.^1 —*■ tt through t-channel stau exchange. 
In addition, the bottom and tau Yukawa couplings increase with tan /3, which increases the annihilation 
cross sections into rs and 5s. Finally, the H and A Higgs boson masses are decreasing with tan/3, and an¬ 
nihilation rates which proceed through these resonances increase. Co-annihilations again gives enhanced 
annihilation cross sections on the left and farthest right hand sides of the allowed parameter space. The 
glitch in contours around TOg ~ 2700 GeV and TOi /2 ~ 425 GeV occurs because to~^ ~ mt = 175 GeV, 

so that a{ZiZi tT) becomes large. In Fig. we show the TOg vs. TOi /2 plane for tan/3 = 45. In this 
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Figure 3: Regions of neutralino relic density in the TOg vs. mi /2 plane for Ag = 0 and tan/3 = 30. 


case, the structure of the plane is changing qualitatively, especially for fj, < 0. First, there is a new region 
of disallowed parameter space for /i < 0 in the lower left due to m\ < 0, which signals a breakdown 
of the REWSB mechanism. Second, a corridor of very low relic density passes diagonally through the 
plot. The center of this region is where 2 to— ~ niA and mn. At the A and H resonance, there is very 
efficient neutralino annihilation into 55 final states. This is illustrated in Fig. H(left), where we show the 
integrated annihilation cross section times velocity versus TOg for TOi /2 = 600 GeV and fi < 0. At the very 
lowest values of TOg, there is again the sharp peak due to neutralino-stau and stau-stau co-annihilations. 
For larger values of TOg, however, the annihilation rate is dominantly into 55 final states over almost the 
entire TOg range. This is due to the large annihilation rates through the s-channel A and H diagrams. 
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Figure 4: Regions of neutralino relic density in th mo vs. mij^ plane for Ap = 0 and tan/3 = 45. 
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Figure 5: Thermally averaged cross section times velocity evaluated at for various subprocesses. The 
thick light-grey(light-blue) curve denotes the total of all annihilation and co-annihilation reactions. Left: 
” 1 - 1/2 = 600 GeV, fj, < 0, Aq = 0 and tan/3 = 45. Right: mi /2 = 300 GeV, > 0, Ap = 0 and tan/3 = 45. 

even when the reactions occur off resonance. In this case, the widths of the A and H are so large (both 
~ 10 — 40 GeV across the range in rnp shown) that efficient s-channel annihilation can occur throughout 
considerable part of the parameter space, even when the resonance condition is not exactly fulfilled. The 
resonance annihilation is explicitly displayed in this plot as the annihilation bump at mp just below 1300 
GeV. Another annihilation possibility is that ZiZi bb via t and u channel graphs. In fact, these 
annihilation graphs are enhanced due to the large b Yukawa coupling and decreasing value of , but we 
have checked that the s-channel annihilation is still far the dominant channel. Annihilation into rf is the 
next most likely channel, but is always below the level of annihilation into bb for the parameters shown 
in Fig. ||(left). At even higher values of mo where the higgsino component of Z\ becomes non-negligible. 
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Figure 6: Neutralino relic density n~ /i^(left) and the fine tuning parameter(right) [defined by Eq.(^)], 
versus ttiq for Aq = 0, TOi /2 = 600 GeV, ^ < 0 and tan/3 = 10, 30 and 45. 



Figure 7: The same as Fig. ^but for mi /2 = 600 GeV and /r < 0. 

the annihilations into WW and ZZ again become important; finally, at the highest values of mo, the Wi 
and Z 2 co-annihilation channels become large. 

In Fig. I^(right), we show again the subprocess annihilation rates versus mo for tan/3 = 45, but this time 
for /r > 0 and for mi /2 = 300 GeV. Although no explicit resonance is evident for /r > 0, the dominant 
annihilations are once again into bb final states over most of the parameter space, due to the wide Higgs 
resonances. To summarize the regions of mSUGRA model parameter space with reasonable values of 
neutralino relic density, we can label four important regions: i.) annihilation through t-channel slepton- 
especially stau- exchange, as occurs for low values of mg and mi/ 2 , ii.) the stau co-annihilation region 
for low values of mg on the edge of the excluded region, iii.) the large mg region with non-negligible 
higgsino-component annihilation, and also Wi (and possibly Z 2 ) co-annihilation occurs near the edge of 
the limit of parameter space, and iv.) annihilation into bb and rr final states through s-channel A and H 
resonances at high tan f3. Other regions can include top or bottom squark co-annihilation for large values 
of Ag, again on the edge of parameter space where ti or bi become light, or annihilation through Z or h 
resonances. The .Z—resonance region is essentially excluded now by constraints on sparticle masses from 
LEP2. It is useful to view the relic density h? directly as a function of model parameters. We show 
in Fig. ||(left) the value of versus the parameter mg for fixed mi/2 = 600 GeV, Ag = 0, /x < 0 

and for tan/3 = 10, 30 and 45. The dashed curves show the result with no co-annihilations, while the 
solid curves yield the complete calculation. The shaded band denotes the cosmologically favored region 
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with 0.1 < r2~ < 0.3. For this value of 1111 / 2 ^ the lower tan/3 curves yield a favored relic density only 


in the very low and very high toq regions, and here the curves have a very sharp slope. The large slope 
is indicative of large fine-tuning, in that a small change of model parameters, in this case toq: yields a 
large change in h^. In contrast, the tan/3 = 45 curve shows a large region with good relic density and 
nearly zero slope(/i > 0) or not very steep slope(/r < 0), and hence with little fine-tuning. In Fig. Il(right), 
we show the corresponding values of the fine-tuning, basically the logarithmic derivative, as advocated 
by Ellis and Olive 


A(mo) = 


Too 


dVt~ h?' 




dmo 


( 8 ) 


As indicated earlier, the low fine-tuning regions mostly coincide with that of neutralino annihilation via 
t-channel slepton exchange (region /.)), or off-resonance annihilation through A and H (region iv.)). 
The co-annihilation region ii.) and focus point region in.) tend to have higher fine-tunings due to the 
steep rise of the cross sections. Regions with simultaneous low fine- tuning and preferred h? values 
are the best candidates for viable mSUGRA parameters. In Fig. ^left), we show Q.~ h? versus toq for 
’ 7 ^l /2 = 300 GeV, Aq = 0, /i > 0 and the same three tan/3 parameters. The curves reflect the broad 
regions of parameter space with reasonable relic density values at high tan /3. The corresponding plot 
of the fine-tuning parameter is shown in Fig. |^(right). Again, there is large fine-tuning at the edges of 
parameter space, but low fine-tuning in the intermediate regions. In conclusion, the relic density and the 
fine-tuning parameter together tend to prefer mSUGRA model parameters in regions i.) or iv.). These 
two regions lead to distinct collider signatures for future searches for supersymmetric matter. 


4 Conclusions 


In conclusion, we have performed a calculation of the neutralino relic density in the minimal supergravity 
model including ^ 12^2 neutralino annihilation and co-annihilation processes, where the initial state 
includes Zi, Z 2 , Wi, ei, /ii, ri, ti and 61 . The calculation was performed using the GompHEP program 
for automatic evaluation of Feynman diagrams, coupled with ISAJET for sparticle mass evaluation in the 
mSUGRA model, and for standard and supersymmetric couplings and decay widths. We implemented 
relativistic thermal averaging, which is especially important for evaluating the relic density when reso¬ 
nances in the annihilation cross section are present, and neutralino thermal velocities can be relativistic. 
The three-dimensional integration was performed by Monte Carlo evaluation with importance sampling, 
which yields in general good convergence even in the presence of narrow resonances. We note that a 
calculation of similar scope and procedure was recently reported in Ref. j^. 

We found four regions of parameter space that led to relic densities in accord with results from cos¬ 
mological measurements, i.e. 0.1 < < 0.3. These include /.) the region dominated by t-channel 

slepton exchange, ii.) the region dominated by stau co-annihilation, Hi.) the large tuq region dominated 
by a more higgsino-like neutralino and iv.) the broad regions at high tan/3 dominated by off-shell annihila¬ 
tion through the A and H Higgs boson resonances. Regions ii.) and Hi.) generally have large fine-tuning 
associated with them, and although it is logically possible that nature has chosen such parameters, any 
slight deviation of model parameters would lead to either too low or too high a relic density. Region 
i.) generally has the property that some of the sleptons have masses less than about 300-400 GeV. This 
region can give rise to a rich set of collider signatures, since many of the sparticles are relatively light. 

Region iv.) gives broad regions of model parameter space with reasonable values of relic density as well 
as low values of the fine-tuning parameter. It can also allow quite heavy values of SUSY particle masses, 
which would be useful to suppress many flavor-violating (such as 6 ^ S 7 )Q and CP violating loop 
processes, and the muon g — 2 value |^. In many respects region iv.) is a favored region of parameter 
space. The neutralino relic density may well point the way to the sort of SUSY signatures we should 
expect at high energy collider experiments. 
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